%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% This program replicates Table 4 in the paper "Optimal
%%% Portfolio Choice with Fat Tails and Parameter Uncertainty",
%%% Journal of Financial and Quantitative Analysis (2024),
%%% Raymond Kan & Nathan Lassance 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

clear
close all

%%% Load data
load Dataset10MOM.txt;
load Dataset16ANOM.txt;
load Dataset25SBETA.txt;
load Dataset25SBTM.txt;
load Dataset25OPINV.txt;
load Dataset30IND.txt

%%% Define variables
gam = 1;
cbps = 10;

%%%%%%%%%%%%%%%%%%
%%% 10MOM dataset
%%%%%%%%%%%%%%%%%%

R = Dataset10MOM;
N = length(R(1,:));
onevec = ones(N,1);
Ttotal = length(R(:,1));

%%% T = 60
T = 60;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T60 = zeros(NumberWindow,N);
wvec2T60 = zeros(NumberWindow,N);
wvec3T60 = zeros(NumberWindow,N);
wvec4T60 = zeros(NumberWindow,N);
Returnvec1T60gross = zeros(NumberWindow,1);
Returnvec2T60gross = zeros(NumberWindow,1);
Returnvec3T60gross = zeros(NumberWindow,1);
Returnvec4T60gross = zeros(NumberWindow,1);
Returnvec1T60net = zeros(NumberWindow,1);
Returnvec2T60net = zeros(NumberWindow,1);
Returnvec3T60net = zeros(NumberWindow,1);
Returnvec4T60net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T60(t,:) = onevec./N;
    Returnvec1T60gross(t,1) = wvec1T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T60(t,:)'-wplus,1);
    end 
    Returnvec1T60net(t,1) = (1+Returnvec1T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T60(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T60gross(t,1) = wvec2T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T60(t,:)'-wplus,1);
    end 
    Returnvec2T60net(t,1) = (1+Returnvec2T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T60(t,:) = wghat;
    Returnvec3T60gross(t,1) = wvec3T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T60(t,:)'-wplus,1);
    end 
    Returnvec3T60net(t,1) = (1+Returnvec3T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T60(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T60gross(t,1) = wvec4T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T60(t,:)'-wplus,1);
    end 
    Returnvec4T60net(t,1) = (1+Returnvec4T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T60gross = 12*(mean(Returnvec1T60gross)-(gam/2)*var(Returnvec1T60gross));
EU2T60gross = 12*(mean(Returnvec2T60gross)-(gam/2)*var(Returnvec2T60gross));
EU3T60gross = 12*(mean(Returnvec3T60gross)-(gam/2)*var(Returnvec3T60gross));
EU4T60gross = 12*(mean(Returnvec4T60gross)-(gam/2)*var(Returnvec4T60gross));
EU1T60net = 12*(mean(Returnvec1T60net)-(gam/2)*var(Returnvec1T60net));
EU2T60net = 12*(mean(Returnvec2T60net)-(gam/2)*var(Returnvec2T60net));
EU3T60net = 12*(mean(Returnvec3T60net)-(gam/2)*var(Returnvec3T60net));
EU4T60net = 12*(mean(Returnvec4T60net)-(gam/2)*var(Returnvec4T60net));

%%% T = 120
T = 120;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T120 = zeros(NumberWindow,N);
wvec2T120 = zeros(NumberWindow,N);
wvec3T120 = zeros(NumberWindow,N);
wvec4T120 = zeros(NumberWindow,N);
Returnvec1T120gross = zeros(NumberWindow,1);
Returnvec2T120gross = zeros(NumberWindow,1);
Returnvec3T120gross = zeros(NumberWindow,1);
Returnvec4T120gross = zeros(NumberWindow,1);
Returnvec1T120net = zeros(NumberWindow,1);
Returnvec2T120net = zeros(NumberWindow,1);
Returnvec3T120net = zeros(NumberWindow,1);
Returnvec4T120net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T120(t,:) = onevec./N;
    Returnvec1T120gross(t,1) = wvec1T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T120(t,:)'-wplus,1);
    end 
    Returnvec1T120net(t,1) = (1+Returnvec1T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T120(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T120gross(t,1) = wvec2T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T120(t,:)'-wplus,1);
    end 
    Returnvec2T120net(t,1) = (1+Returnvec2T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T120(t,:) = wghat;
    Returnvec3T120gross(t,1) = wvec3T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T120(t,:)'-wplus,1);
    end 
    Returnvec3T120net(t,1) = (1+Returnvec3T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T120(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T120gross(t,1) = wvec4T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T120(t,:)'-wplus,1);
    end 
    Returnvec4T120net(t,1) = (1+Returnvec4T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T120gross = 12*(mean(Returnvec1T120gross)-(gam/2)*var(Returnvec1T120gross));
EU2T120gross = 12*(mean(Returnvec2T120gross)-(gam/2)*var(Returnvec2T120gross));
EU3T120gross = 12*(mean(Returnvec3T120gross)-(gam/2)*var(Returnvec3T120gross));
EU4T120gross = 12*(mean(Returnvec4T120gross)-(gam/2)*var(Returnvec4T120gross));
EU1T120net = 12*(mean(Returnvec1T120net)-(gam/2)*var(Returnvec1T120net));
EU2T120net = 12*(mean(Returnvec2T120net)-(gam/2)*var(Returnvec2T120net));
EU3T120net = 12*(mean(Returnvec3T120net)-(gam/2)*var(Returnvec3T120net));
EU4T120net = 12*(mean(Returnvec4T120net)-(gam/2)*var(Returnvec4T120net));

%%% T = 240
T = 240;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T240 = zeros(NumberWindow,N);
wvec2T240 = zeros(NumberWindow,N);
wvec3T240 = zeros(NumberWindow,N);
wvec4T240 = zeros(NumberWindow,N);
Returnvec1T240gross = zeros(NumberWindow,1);
Returnvec2T240gross = zeros(NumberWindow,1);
Returnvec3T240gross = zeros(NumberWindow,1);
Returnvec4T240gross = zeros(NumberWindow,1);
Returnvec1T240net = zeros(NumberWindow,1);
Returnvec2T240net = zeros(NumberWindow,1);
Returnvec3T240net = zeros(NumberWindow,1);
Returnvec4T240net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T240(t,:) = onevec./N;
    Returnvec1T240gross(t,1) = wvec1T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T240(t,:)'-wplus,1);
    end 
    Returnvec1T240net(t,1) = (1+Returnvec1T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T240(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T240gross(t,1) = wvec2T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T240(t,:)'-wplus,1);
    end 
    Returnvec2T240net(t,1) = (1+Returnvec2T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T240(t,:) = wghat;
    Returnvec3T240gross(t,1) = wvec3T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T240(t,:)'-wplus,1);
    end 
    Returnvec3T240net(t,1) = (1+Returnvec3T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T240(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T240gross(t,1) = wvec4T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T240(t,:)'-wplus,1);
    end 
    Returnvec4T240net(t,1) = (1+Returnvec4T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T240gross = 12*(mean(Returnvec1T240gross)-(gam/2)*var(Returnvec1T240gross));
EU2T240gross = 12*(mean(Returnvec2T240gross)-(gam/2)*var(Returnvec2T240gross));
EU3T240gross = 12*(mean(Returnvec3T240gross)-(gam/2)*var(Returnvec3T240gross));
EU4T240gross = 12*(mean(Returnvec4T240gross)-(gam/2)*var(Returnvec4T240gross));
EU1T240net = 12*(mean(Returnvec1T240net)-(gam/2)*var(Returnvec1T240net));
EU2T240net = 12*(mean(Returnvec2T240net)-(gam/2)*var(Returnvec2T240net));
EU3T240net = 12*(mean(Returnvec3T240net)-(gam/2)*var(Returnvec3T240net));
EU4T240net = 12*(mean(Returnvec4T240net)-(gam/2)*var(Returnvec4T240net));

Table4Panel10MOM = ...
[EU1T60gross,EU1T60net,EU1T120gross,EU1T120net,EU1T240gross,EU1T240net;...
EU2T60gross,EU2T60net,EU2T120gross,EU2T120net,EU2T240gross,EU2T240net;...
EU3T60gross,EU3T60net,EU3T120gross,EU3T120net,EU3T240gross,EU3T240net;...
EU4T60gross,EU4T60net,EU4T120gross,EU4T120net,EU4T240gross,EU4T240net]


%%%%%%%%%%%%%%%%%%
%%% 16ANOM dataset
%%%%%%%%%%%%%%%%%%

R = Dataset16ANOM;
N = length(R(1,:));
onevec = ones(N,1);
Ttotal = length(R(:,1));

%%% T = 60
T = 60;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T60 = zeros(NumberWindow,N);
wvec2T60 = zeros(NumberWindow,N);
wvec3T60 = zeros(NumberWindow,N);
wvec4T60 = zeros(NumberWindow,N);
Returnvec1T60gross = zeros(NumberWindow,1);
Returnvec2T60gross = zeros(NumberWindow,1);
Returnvec3T60gross = zeros(NumberWindow,1);
Returnvec4T60gross = zeros(NumberWindow,1);
Returnvec1T60net = zeros(NumberWindow,1);
Returnvec2T60net = zeros(NumberWindow,1);
Returnvec3T60net = zeros(NumberWindow,1);
Returnvec4T60net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T60(t,:) = onevec./N;
    Returnvec1T60gross(t,1) = wvec1T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T60(t,:)'-wplus,1);
    end 
    Returnvec1T60net(t,1) = (1+Returnvec1T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T60(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T60gross(t,1) = wvec2T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T60(t,:)'-wplus,1);
    end 
    Returnvec2T60net(t,1) = (1+Returnvec2T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T60(t,:) = wghat;
    Returnvec3T60gross(t,1) = wvec3T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T60(t,:)'-wplus,1);
    end 
    Returnvec3T60net(t,1) = (1+Returnvec3T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T60(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T60gross(t,1) = wvec4T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T60(t,:)'-wplus,1);
    end 
    Returnvec4T60net(t,1) = (1+Returnvec4T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T60gross = 12*(mean(Returnvec1T60gross)-(gam/2)*var(Returnvec1T60gross));
EU2T60gross = 12*(mean(Returnvec2T60gross)-(gam/2)*var(Returnvec2T60gross));
EU3T60gross = 12*(mean(Returnvec3T60gross)-(gam/2)*var(Returnvec3T60gross));
EU4T60gross = 12*(mean(Returnvec4T60gross)-(gam/2)*var(Returnvec4T60gross));
EU1T60net = 12*(mean(Returnvec1T60net)-(gam/2)*var(Returnvec1T60net));
EU2T60net = 12*(mean(Returnvec2T60net)-(gam/2)*var(Returnvec2T60net));
EU3T60net = 12*(mean(Returnvec3T60net)-(gam/2)*var(Returnvec3T60net));
EU4T60net = 12*(mean(Returnvec4T60net)-(gam/2)*var(Returnvec4T60net));

%%% T = 120
T = 120;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T120 = zeros(NumberWindow,N);
wvec2T120 = zeros(NumberWindow,N);
wvec3T120 = zeros(NumberWindow,N);
wvec4T120 = zeros(NumberWindow,N);
Returnvec1T120gross = zeros(NumberWindow,1);
Returnvec2T120gross = zeros(NumberWindow,1);
Returnvec3T120gross = zeros(NumberWindow,1);
Returnvec4T120gross = zeros(NumberWindow,1);
Returnvec1T120net = zeros(NumberWindow,1);
Returnvec2T120net = zeros(NumberWindow,1);
Returnvec3T120net = zeros(NumberWindow,1);
Returnvec4T120net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T120(t,:) = onevec./N;
    Returnvec1T120gross(t,1) = wvec1T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T120(t,:)'-wplus,1);
    end 
    Returnvec1T120net(t,1) = (1+Returnvec1T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T120(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T120gross(t,1) = wvec2T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T120(t,:)'-wplus,1);
    end 
    Returnvec2T120net(t,1) = (1+Returnvec2T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T120(t,:) = wghat;
    Returnvec3T120gross(t,1) = wvec3T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T120(t,:)'-wplus,1);
    end 
    Returnvec3T120net(t,1) = (1+Returnvec3T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T120(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T120gross(t,1) = wvec4T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T120(t,:)'-wplus,1);
    end 
    Returnvec4T120net(t,1) = (1+Returnvec4T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T120gross = 12*(mean(Returnvec1T120gross)-(gam/2)*var(Returnvec1T120gross));
EU2T120gross = 12*(mean(Returnvec2T120gross)-(gam/2)*var(Returnvec2T120gross));
EU3T120gross = 12*(mean(Returnvec3T120gross)-(gam/2)*var(Returnvec3T120gross));
EU4T120gross = 12*(mean(Returnvec4T120gross)-(gam/2)*var(Returnvec4T120gross));
EU1T120net = 12*(mean(Returnvec1T120net)-(gam/2)*var(Returnvec1T120net));
EU2T120net = 12*(mean(Returnvec2T120net)-(gam/2)*var(Returnvec2T120net));
EU3T120net = 12*(mean(Returnvec3T120net)-(gam/2)*var(Returnvec3T120net));
EU4T120net = 12*(mean(Returnvec4T120net)-(gam/2)*var(Returnvec4T120net));

%%% T = 240
T = 240;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T240 = zeros(NumberWindow,N);
wvec2T240 = zeros(NumberWindow,N);
wvec3T240 = zeros(NumberWindow,N);
wvec4T240 = zeros(NumberWindow,N);
Returnvec1T240gross = zeros(NumberWindow,1);
Returnvec2T240gross = zeros(NumberWindow,1);
Returnvec3T240gross = zeros(NumberWindow,1);
Returnvec4T240gross = zeros(NumberWindow,1);
Returnvec1T240net = zeros(NumberWindow,1);
Returnvec2T240net = zeros(NumberWindow,1);
Returnvec3T240net = zeros(NumberWindow,1);
Returnvec4T240net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T240(t,:) = onevec./N;
    Returnvec1T240gross(t,1) = wvec1T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T240(t,:)'-wplus,1);
    end 
    Returnvec1T240net(t,1) = (1+Returnvec1T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T240(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T240gross(t,1) = wvec2T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T240(t,:)'-wplus,1);
    end 
    Returnvec2T240net(t,1) = (1+Returnvec2T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T240(t,:) = wghat;
    Returnvec3T240gross(t,1) = wvec3T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T240(t,:)'-wplus,1);
    end 
    Returnvec3T240net(t,1) = (1+Returnvec3T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T240(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T240gross(t,1) = wvec4T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T240(t,:)'-wplus,1);
    end 
    Returnvec4T240net(t,1) = (1+Returnvec4T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T240gross = 12*(mean(Returnvec1T240gross)-(gam/2)*var(Returnvec1T240gross));
EU2T240gross = 12*(mean(Returnvec2T240gross)-(gam/2)*var(Returnvec2T240gross));
EU3T240gross = 12*(mean(Returnvec3T240gross)-(gam/2)*var(Returnvec3T240gross));
EU4T240gross = 12*(mean(Returnvec4T240gross)-(gam/2)*var(Returnvec4T240gross));
EU1T240net = 12*(mean(Returnvec1T240net)-(gam/2)*var(Returnvec1T240net));
EU2T240net = 12*(mean(Returnvec2T240net)-(gam/2)*var(Returnvec2T240net));
EU3T240net = 12*(mean(Returnvec3T240net)-(gam/2)*var(Returnvec3T240net));
EU4T240net = 12*(mean(Returnvec4T240net)-(gam/2)*var(Returnvec4T240net));

Table4Panel16ANOM = ...
[EU1T60gross,EU1T60net,EU1T120gross,EU1T120net,EU1T240gross,EU1T240net;...
EU2T60gross,EU2T60net,EU2T120gross,EU2T120net,EU2T240gross,EU2T240net;...
EU3T60gross,EU3T60net,EU3T120gross,EU3T120net,EU3T240gross,EU3T240net;...
EU4T60gross,EU4T60net,EU4T120gross,EU4T120net,EU4T240gross,EU4T240net]


%%%%%%%%%%%%%%%%%%
%%% 25SBETA dataset
%%%%%%%%%%%%%%%%%%

R = Dataset25SBETA;
N = length(R(1,:));
onevec = ones(N,1);
Ttotal = length(R(:,1));

%%% T = 60
T = 60;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T60 = zeros(NumberWindow,N);
wvec2T60 = zeros(NumberWindow,N);
wvec3T60 = zeros(NumberWindow,N);
wvec4T60 = zeros(NumberWindow,N);
Returnvec1T60gross = zeros(NumberWindow,1);
Returnvec2T60gross = zeros(NumberWindow,1);
Returnvec3T60gross = zeros(NumberWindow,1);
Returnvec4T60gross = zeros(NumberWindow,1);
Returnvec1T60net = zeros(NumberWindow,1);
Returnvec2T60net = zeros(NumberWindow,1);
Returnvec3T60net = zeros(NumberWindow,1);
Returnvec4T60net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T60(t,:) = onevec./N;
    Returnvec1T60gross(t,1) = wvec1T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T60(t,:)'-wplus,1);
    end 
    Returnvec1T60net(t,1) = (1+Returnvec1T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T60(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T60gross(t,1) = wvec2T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T60(t,:)'-wplus,1);
    end 
    Returnvec2T60net(t,1) = (1+Returnvec2T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T60(t,:) = wghat;
    Returnvec3T60gross(t,1) = wvec3T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T60(t,:)'-wplus,1);
    end 
    Returnvec3T60net(t,1) = (1+Returnvec3T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T60(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T60gross(t,1) = wvec4T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T60(t,:)'-wplus,1);
    end 
    Returnvec4T60net(t,1) = (1+Returnvec4T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T60gross = 12*(mean(Returnvec1T60gross)-(gam/2)*var(Returnvec1T60gross));
EU2T60gross = 12*(mean(Returnvec2T60gross)-(gam/2)*var(Returnvec2T60gross));
EU3T60gross = 12*(mean(Returnvec3T60gross)-(gam/2)*var(Returnvec3T60gross));
EU4T60gross = 12*(mean(Returnvec4T60gross)-(gam/2)*var(Returnvec4T60gross));
EU1T60net = 12*(mean(Returnvec1T60net)-(gam/2)*var(Returnvec1T60net));
EU2T60net = 12*(mean(Returnvec2T60net)-(gam/2)*var(Returnvec2T60net));
EU3T60net = 12*(mean(Returnvec3T60net)-(gam/2)*var(Returnvec3T60net));
EU4T60net = 12*(mean(Returnvec4T60net)-(gam/2)*var(Returnvec4T60net));

%%% T = 120
T = 120;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T120 = zeros(NumberWindow,N);
wvec2T120 = zeros(NumberWindow,N);
wvec3T120 = zeros(NumberWindow,N);
wvec4T120 = zeros(NumberWindow,N);
Returnvec1T120gross = zeros(NumberWindow,1);
Returnvec2T120gross = zeros(NumberWindow,1);
Returnvec3T120gross = zeros(NumberWindow,1);
Returnvec4T120gross = zeros(NumberWindow,1);
Returnvec1T120net = zeros(NumberWindow,1);
Returnvec2T120net = zeros(NumberWindow,1);
Returnvec3T120net = zeros(NumberWindow,1);
Returnvec4T120net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T120(t,:) = onevec./N;
    Returnvec1T120gross(t,1) = wvec1T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T120(t,:)'-wplus,1);
    end 
    Returnvec1T120net(t,1) = (1+Returnvec1T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T120(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T120gross(t,1) = wvec2T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T120(t,:)'-wplus,1);
    end 
    Returnvec2T120net(t,1) = (1+Returnvec2T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T120(t,:) = wghat;
    Returnvec3T120gross(t,1) = wvec3T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T120(t,:)'-wplus,1);
    end 
    Returnvec3T120net(t,1) = (1+Returnvec3T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T120(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T120gross(t,1) = wvec4T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T120(t,:)'-wplus,1);
    end 
    Returnvec4T120net(t,1) = (1+Returnvec4T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T120gross = 12*(mean(Returnvec1T120gross)-(gam/2)*var(Returnvec1T120gross));
EU2T120gross = 12*(mean(Returnvec2T120gross)-(gam/2)*var(Returnvec2T120gross));
EU3T120gross = 12*(mean(Returnvec3T120gross)-(gam/2)*var(Returnvec3T120gross));
EU4T120gross = 12*(mean(Returnvec4T120gross)-(gam/2)*var(Returnvec4T120gross));
EU1T120net = 12*(mean(Returnvec1T120net)-(gam/2)*var(Returnvec1T120net));
EU2T120net = 12*(mean(Returnvec2T120net)-(gam/2)*var(Returnvec2T120net));
EU3T120net = 12*(mean(Returnvec3T120net)-(gam/2)*var(Returnvec3T120net));
EU4T120net = 12*(mean(Returnvec4T120net)-(gam/2)*var(Returnvec4T120net));

%%% T = 240
T = 240;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T240 = zeros(NumberWindow,N);
wvec2T240 = zeros(NumberWindow,N);
wvec3T240 = zeros(NumberWindow,N);
wvec4T240 = zeros(NumberWindow,N);
Returnvec1T240gross = zeros(NumberWindow,1);
Returnvec2T240gross = zeros(NumberWindow,1);
Returnvec3T240gross = zeros(NumberWindow,1);
Returnvec4T240gross = zeros(NumberWindow,1);
Returnvec1T240net = zeros(NumberWindow,1);
Returnvec2T240net = zeros(NumberWindow,1);
Returnvec3T240net = zeros(NumberWindow,1);
Returnvec4T240net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T240(t,:) = onevec./N;
    Returnvec1T240gross(t,1) = wvec1T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T240(t,:)'-wplus,1);
    end 
    Returnvec1T240net(t,1) = (1+Returnvec1T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T240(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T240gross(t,1) = wvec2T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T240(t,:)'-wplus,1);
    end 
    Returnvec2T240net(t,1) = (1+Returnvec2T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T240(t,:) = wghat;
    Returnvec3T240gross(t,1) = wvec3T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T240(t,:)'-wplus,1);
    end 
    Returnvec3T240net(t,1) = (1+Returnvec3T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T240(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T240gross(t,1) = wvec4T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T240(t,:)'-wplus,1);
    end 
    Returnvec4T240net(t,1) = (1+Returnvec4T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T240gross = 12*(mean(Returnvec1T240gross)-(gam/2)*var(Returnvec1T240gross));
EU2T240gross = 12*(mean(Returnvec2T240gross)-(gam/2)*var(Returnvec2T240gross));
EU3T240gross = 12*(mean(Returnvec3T240gross)-(gam/2)*var(Returnvec3T240gross));
EU4T240gross = 12*(mean(Returnvec4T240gross)-(gam/2)*var(Returnvec4T240gross));
EU1T240net = 12*(mean(Returnvec1T240net)-(gam/2)*var(Returnvec1T240net));
EU2T240net = 12*(mean(Returnvec2T240net)-(gam/2)*var(Returnvec2T240net));
EU3T240net = 12*(mean(Returnvec3T240net)-(gam/2)*var(Returnvec3T240net));
EU4T240net = 12*(mean(Returnvec4T240net)-(gam/2)*var(Returnvec4T240net));

Table4Panel25SBETA = ...
[EU1T60gross,EU1T60net,EU1T120gross,EU1T120net,EU1T240gross,EU1T240net;...
EU2T60gross,EU2T60net,EU2T120gross,EU2T120net,EU2T240gross,EU2T240net;...
EU3T60gross,EU3T60net,EU3T120gross,EU3T120net,EU3T240gross,EU3T240net;...
EU4T60gross,EU4T60net,EU4T120gross,EU4T120net,EU4T240gross,EU4T240net]


%%%%%%%%%%%%%%%%%%
%%% 25SBTM dataset
%%%%%%%%%%%%%%%%%%

R = Dataset25SBTM;
N = length(R(1,:));
onevec = ones(N,1);
Ttotal = length(R(:,1));

%%% T = 60
T = 60;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T60 = zeros(NumberWindow,N);
wvec2T60 = zeros(NumberWindow,N);
wvec3T60 = zeros(NumberWindow,N);
wvec4T60 = zeros(NumberWindow,N);
Returnvec1T60gross = zeros(NumberWindow,1);
Returnvec2T60gross = zeros(NumberWindow,1);
Returnvec3T60gross = zeros(NumberWindow,1);
Returnvec4T60gross = zeros(NumberWindow,1);
Returnvec1T60net = zeros(NumberWindow,1);
Returnvec2T60net = zeros(NumberWindow,1);
Returnvec3T60net = zeros(NumberWindow,1);
Returnvec4T60net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T60(t,:) = onevec./N;
    Returnvec1T60gross(t,1) = wvec1T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T60(t,:)'-wplus,1);
    end 
    Returnvec1T60net(t,1) = (1+Returnvec1T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T60(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T60gross(t,1) = wvec2T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T60(t,:)'-wplus,1);
    end 
    Returnvec2T60net(t,1) = (1+Returnvec2T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T60(t,:) = wghat;
    Returnvec3T60gross(t,1) = wvec3T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T60(t,:)'-wplus,1);
    end 
    Returnvec3T60net(t,1) = (1+Returnvec3T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T60(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T60gross(t,1) = wvec4T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T60(t,:)'-wplus,1);
    end 
    Returnvec4T60net(t,1) = (1+Returnvec4T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T60gross = 12*(mean(Returnvec1T60gross)-(gam/2)*var(Returnvec1T60gross));
EU2T60gross = 12*(mean(Returnvec2T60gross)-(gam/2)*var(Returnvec2T60gross));
EU3T60gross = 12*(mean(Returnvec3T60gross)-(gam/2)*var(Returnvec3T60gross));
EU4T60gross = 12*(mean(Returnvec4T60gross)-(gam/2)*var(Returnvec4T60gross));
EU1T60net = 12*(mean(Returnvec1T60net)-(gam/2)*var(Returnvec1T60net));
EU2T60net = 12*(mean(Returnvec2T60net)-(gam/2)*var(Returnvec2T60net));
EU3T60net = 12*(mean(Returnvec3T60net)-(gam/2)*var(Returnvec3T60net));
EU4T60net = 12*(mean(Returnvec4T60net)-(gam/2)*var(Returnvec4T60net));

%%% T = 120
T = 120;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T120 = zeros(NumberWindow,N);
wvec2T120 = zeros(NumberWindow,N);
wvec3T120 = zeros(NumberWindow,N);
wvec4T120 = zeros(NumberWindow,N);
Returnvec1T120gross = zeros(NumberWindow,1);
Returnvec2T120gross = zeros(NumberWindow,1);
Returnvec3T120gross = zeros(NumberWindow,1);
Returnvec4T120gross = zeros(NumberWindow,1);
Returnvec1T120net = zeros(NumberWindow,1);
Returnvec2T120net = zeros(NumberWindow,1);
Returnvec3T120net = zeros(NumberWindow,1);
Returnvec4T120net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T120(t,:) = onevec./N;
    Returnvec1T120gross(t,1) = wvec1T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T120(t,:)'-wplus,1);
    end 
    Returnvec1T120net(t,1) = (1+Returnvec1T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T120(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T120gross(t,1) = wvec2T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T120(t,:)'-wplus,1);
    end 
    Returnvec2T120net(t,1) = (1+Returnvec2T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T120(t,:) = wghat;
    Returnvec3T120gross(t,1) = wvec3T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T120(t,:)'-wplus,1);
    end 
    Returnvec3T120net(t,1) = (1+Returnvec3T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T120(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T120gross(t,1) = wvec4T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T120(t,:)'-wplus,1);
    end 
    Returnvec4T120net(t,1) = (1+Returnvec4T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T120gross = 12*(mean(Returnvec1T120gross)-(gam/2)*var(Returnvec1T120gross));
EU2T120gross = 12*(mean(Returnvec2T120gross)-(gam/2)*var(Returnvec2T120gross));
EU3T120gross = 12*(mean(Returnvec3T120gross)-(gam/2)*var(Returnvec3T120gross));
EU4T120gross = 12*(mean(Returnvec4T120gross)-(gam/2)*var(Returnvec4T120gross));
EU1T120net = 12*(mean(Returnvec1T120net)-(gam/2)*var(Returnvec1T120net));
EU2T120net = 12*(mean(Returnvec2T120net)-(gam/2)*var(Returnvec2T120net));
EU3T120net = 12*(mean(Returnvec3T120net)-(gam/2)*var(Returnvec3T120net));
EU4T120net = 12*(mean(Returnvec4T120net)-(gam/2)*var(Returnvec4T120net));

%%% T = 240
T = 240;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T240 = zeros(NumberWindow,N);
wvec2T240 = zeros(NumberWindow,N);
wvec3T240 = zeros(NumberWindow,N);
wvec4T240 = zeros(NumberWindow,N);
Returnvec1T240gross = zeros(NumberWindow,1);
Returnvec2T240gross = zeros(NumberWindow,1);
Returnvec3T240gross = zeros(NumberWindow,1);
Returnvec4T240gross = zeros(NumberWindow,1);
Returnvec1T240net = zeros(NumberWindow,1);
Returnvec2T240net = zeros(NumberWindow,1);
Returnvec3T240net = zeros(NumberWindow,1);
Returnvec4T240net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T240(t,:) = onevec./N;
    Returnvec1T240gross(t,1) = wvec1T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T240(t,:)'-wplus,1);
    end 
    Returnvec1T240net(t,1) = (1+Returnvec1T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T240(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T240gross(t,1) = wvec2T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T240(t,:)'-wplus,1);
    end 
    Returnvec2T240net(t,1) = (1+Returnvec2T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T240(t,:) = wghat;
    Returnvec3T240gross(t,1) = wvec3T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T240(t,:)'-wplus,1);
    end 
    Returnvec3T240net(t,1) = (1+Returnvec3T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T240(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T240gross(t,1) = wvec4T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T240(t,:)'-wplus,1);
    end 
    Returnvec4T240net(t,1) = (1+Returnvec4T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T240gross = 12*(mean(Returnvec1T240gross)-(gam/2)*var(Returnvec1T240gross));
EU2T240gross = 12*(mean(Returnvec2T240gross)-(gam/2)*var(Returnvec2T240gross));
EU3T240gross = 12*(mean(Returnvec3T240gross)-(gam/2)*var(Returnvec3T240gross));
EU4T240gross = 12*(mean(Returnvec4T240gross)-(gam/2)*var(Returnvec4T240gross));
EU1T240net = 12*(mean(Returnvec1T240net)-(gam/2)*var(Returnvec1T240net));
EU2T240net = 12*(mean(Returnvec2T240net)-(gam/2)*var(Returnvec2T240net));
EU3T240net = 12*(mean(Returnvec3T240net)-(gam/2)*var(Returnvec3T240net));
EU4T240net = 12*(mean(Returnvec4T240net)-(gam/2)*var(Returnvec4T240net));

Table4Panel25SBTM = ...
[EU1T60gross,EU1T60net,EU1T120gross,EU1T120net,EU1T240gross,EU1T240net;...
EU2T60gross,EU2T60net,EU2T120gross,EU2T120net,EU2T240gross,EU2T240net;...
EU3T60gross,EU3T60net,EU3T120gross,EU3T120net,EU3T240gross,EU3T240net;...
EU4T60gross,EU4T60net,EU4T120gross,EU4T120net,EU4T240gross,EU4T240net]


%%%%%%%%%%%%%%%%%%
%%% 25OPINV dataset
%%%%%%%%%%%%%%%%%%

R = Dataset25OPINV;
N = length(R(1,:));
onevec = ones(N,1);
Ttotal = length(R(:,1));

%%% T = 60
T = 60;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T60 = zeros(NumberWindow,N);
wvec2T60 = zeros(NumberWindow,N);
wvec3T60 = zeros(NumberWindow,N);
wvec4T60 = zeros(NumberWindow,N);
Returnvec1T60gross = zeros(NumberWindow,1);
Returnvec2T60gross = zeros(NumberWindow,1);
Returnvec3T60gross = zeros(NumberWindow,1);
Returnvec4T60gross = zeros(NumberWindow,1);
Returnvec1T60net = zeros(NumberWindow,1);
Returnvec2T60net = zeros(NumberWindow,1);
Returnvec3T60net = zeros(NumberWindow,1);
Returnvec4T60net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T60(t,:) = onevec./N;
    Returnvec1T60gross(t,1) = wvec1T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T60(t,:)'-wplus,1);
    end 
    Returnvec1T60net(t,1) = (1+Returnvec1T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T60(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T60gross(t,1) = wvec2T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T60(t,:)'-wplus,1);
    end 
    Returnvec2T60net(t,1) = (1+Returnvec2T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T60(t,:) = wghat;
    Returnvec3T60gross(t,1) = wvec3T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T60(t,:)'-wplus,1);
    end 
    Returnvec3T60net(t,1) = (1+Returnvec3T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T60(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T60gross(t,1) = wvec4T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T60(t,:)'-wplus,1);
    end 
    Returnvec4T60net(t,1) = (1+Returnvec4T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T60gross = 12*(mean(Returnvec1T60gross)-(gam/2)*var(Returnvec1T60gross));
EU2T60gross = 12*(mean(Returnvec2T60gross)-(gam/2)*var(Returnvec2T60gross));
EU3T60gross = 12*(mean(Returnvec3T60gross)-(gam/2)*var(Returnvec3T60gross));
EU4T60gross = 12*(mean(Returnvec4T60gross)-(gam/2)*var(Returnvec4T60gross));
EU1T60net = 12*(mean(Returnvec1T60net)-(gam/2)*var(Returnvec1T60net));
EU2T60net = 12*(mean(Returnvec2T60net)-(gam/2)*var(Returnvec2T60net));
EU3T60net = 12*(mean(Returnvec3T60net)-(gam/2)*var(Returnvec3T60net));
EU4T60net = 12*(mean(Returnvec4T60net)-(gam/2)*var(Returnvec4T60net));

%%% T = 120
T = 120;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T120 = zeros(NumberWindow,N);
wvec2T120 = zeros(NumberWindow,N);
wvec3T120 = zeros(NumberWindow,N);
wvec4T120 = zeros(NumberWindow,N);
Returnvec1T120gross = zeros(NumberWindow,1);
Returnvec2T120gross = zeros(NumberWindow,1);
Returnvec3T120gross = zeros(NumberWindow,1);
Returnvec4T120gross = zeros(NumberWindow,1);
Returnvec1T120net = zeros(NumberWindow,1);
Returnvec2T120net = zeros(NumberWindow,1);
Returnvec3T120net = zeros(NumberWindow,1);
Returnvec4T120net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T120(t,:) = onevec./N;
    Returnvec1T120gross(t,1) = wvec1T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T120(t,:)'-wplus,1);
    end 
    Returnvec1T120net(t,1) = (1+Returnvec1T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T120(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T120gross(t,1) = wvec2T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T120(t,:)'-wplus,1);
    end 
    Returnvec2T120net(t,1) = (1+Returnvec2T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T120(t,:) = wghat;
    Returnvec3T120gross(t,1) = wvec3T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T120(t,:)'-wplus,1);
    end 
    Returnvec3T120net(t,1) = (1+Returnvec3T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T120(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T120gross(t,1) = wvec4T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T120(t,:)'-wplus,1);
    end 
    Returnvec4T120net(t,1) = (1+Returnvec4T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T120gross = 12*(mean(Returnvec1T120gross)-(gam/2)*var(Returnvec1T120gross));
EU2T120gross = 12*(mean(Returnvec2T120gross)-(gam/2)*var(Returnvec2T120gross));
EU3T120gross = 12*(mean(Returnvec3T120gross)-(gam/2)*var(Returnvec3T120gross));
EU4T120gross = 12*(mean(Returnvec4T120gross)-(gam/2)*var(Returnvec4T120gross));
EU1T120net = 12*(mean(Returnvec1T120net)-(gam/2)*var(Returnvec1T120net));
EU2T120net = 12*(mean(Returnvec2T120net)-(gam/2)*var(Returnvec2T120net));
EU3T120net = 12*(mean(Returnvec3T120net)-(gam/2)*var(Returnvec3T120net));
EU4T120net = 12*(mean(Returnvec4T120net)-(gam/2)*var(Returnvec4T120net));

%%% T = 240
T = 240;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T240 = zeros(NumberWindow,N);
wvec2T240 = zeros(NumberWindow,N);
wvec3T240 = zeros(NumberWindow,N);
wvec4T240 = zeros(NumberWindow,N);
Returnvec1T240gross = zeros(NumberWindow,1);
Returnvec2T240gross = zeros(NumberWindow,1);
Returnvec3T240gross = zeros(NumberWindow,1);
Returnvec4T240gross = zeros(NumberWindow,1);
Returnvec1T240net = zeros(NumberWindow,1);
Returnvec2T240net = zeros(NumberWindow,1);
Returnvec3T240net = zeros(NumberWindow,1);
Returnvec4T240net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T240(t,:) = onevec./N;
    Returnvec1T240gross(t,1) = wvec1T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T240(t,:)'-wplus,1);
    end 
    Returnvec1T240net(t,1) = (1+Returnvec1T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T240(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T240gross(t,1) = wvec2T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T240(t,:)'-wplus,1);
    end 
    Returnvec2T240net(t,1) = (1+Returnvec2T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T240(t,:) = wghat;
    Returnvec3T240gross(t,1) = wvec3T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T240(t,:)'-wplus,1);
    end 
    Returnvec3T240net(t,1) = (1+Returnvec3T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T240(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T240gross(t,1) = wvec4T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T240(t,:)'-wplus,1);
    end 
    Returnvec4T240net(t,1) = (1+Returnvec4T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T240gross = 12*(mean(Returnvec1T240gross)-(gam/2)*var(Returnvec1T240gross));
EU2T240gross = 12*(mean(Returnvec2T240gross)-(gam/2)*var(Returnvec2T240gross));
EU3T240gross = 12*(mean(Returnvec3T240gross)-(gam/2)*var(Returnvec3T240gross));
EU4T240gross = 12*(mean(Returnvec4T240gross)-(gam/2)*var(Returnvec4T240gross));
EU1T240net = 12*(mean(Returnvec1T240net)-(gam/2)*var(Returnvec1T240net));
EU2T240net = 12*(mean(Returnvec2T240net)-(gam/2)*var(Returnvec2T240net));
EU3T240net = 12*(mean(Returnvec3T240net)-(gam/2)*var(Returnvec3T240net));
EU4T240net = 12*(mean(Returnvec4T240net)-(gam/2)*var(Returnvec4T240net));

Table4Panel25OPINV = ...
[EU1T60gross,EU1T60net,EU1T120gross,EU1T120net,EU1T240gross,EU1T240net;...
EU2T60gross,EU2T60net,EU2T120gross,EU2T120net,EU2T240gross,EU2T240net;...
EU3T60gross,EU3T60net,EU3T120gross,EU3T120net,EU3T240gross,EU3T240net;...
EU4T60gross,EU4T60net,EU4T120gross,EU4T120net,EU4T240gross,EU4T240net]


%%%%%%%%%%%%%%%%%%
%%% 30IND dataset
%%%%%%%%%%%%%%%%%%

R = Dataset30IND;
N = length(R(1,:));
onevec = ones(N,1);
Ttotal = length(R(:,1));

%%% T = 60
T = 60;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T60 = zeros(NumberWindow,N);
wvec2T60 = zeros(NumberWindow,N);
wvec3T60 = zeros(NumberWindow,N);
wvec4T60 = zeros(NumberWindow,N);
Returnvec1T60gross = zeros(NumberWindow,1);
Returnvec2T60gross = zeros(NumberWindow,1);
Returnvec3T60gross = zeros(NumberWindow,1);
Returnvec4T60gross = zeros(NumberWindow,1);
Returnvec1T60net = zeros(NumberWindow,1);
Returnvec2T60net = zeros(NumberWindow,1);
Returnvec3T60net = zeros(NumberWindow,1);
Returnvec4T60net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T60(t,:) = onevec./N;
    Returnvec1T60gross(t,1) = wvec1T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T60(t,:)'-wplus,1);
    end 
    Returnvec1T60net(t,1) = (1+Returnvec1T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T60(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T60gross(t,1) = wvec2T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T60(t,:)'-wplus,1);
    end 
    Returnvec2T60net(t,1) = (1+Returnvec2T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T60(t,:) = wghat;
    Returnvec3T60gross(t,1) = wvec3T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T60(t,:)'-wplus,1);
    end 
    Returnvec3T60net(t,1) = (1+Returnvec3T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T60(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T60gross(t,1) = wvec4T60(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T60(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T60(t,:)'-wplus,1);
    end 
    Returnvec4T60net(t,1) = (1+Returnvec4T60gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T60gross = 12*(mean(Returnvec1T60gross)-(gam/2)*var(Returnvec1T60gross));
EU2T60gross = 12*(mean(Returnvec2T60gross)-(gam/2)*var(Returnvec2T60gross));
EU3T60gross = 12*(mean(Returnvec3T60gross)-(gam/2)*var(Returnvec3T60gross));
EU4T60gross = 12*(mean(Returnvec4T60gross)-(gam/2)*var(Returnvec4T60gross));
EU1T60net = 12*(mean(Returnvec1T60net)-(gam/2)*var(Returnvec1T60net));
EU2T60net = 12*(mean(Returnvec2T60net)-(gam/2)*var(Returnvec2T60net));
EU3T60net = 12*(mean(Returnvec3T60net)-(gam/2)*var(Returnvec3T60net));
EU4T60net = 12*(mean(Returnvec4T60net)-(gam/2)*var(Returnvec4T60net));

%%% T = 120
T = 120;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T120 = zeros(NumberWindow,N);
wvec2T120 = zeros(NumberWindow,N);
wvec3T120 = zeros(NumberWindow,N);
wvec4T120 = zeros(NumberWindow,N);
Returnvec1T120gross = zeros(NumberWindow,1);
Returnvec2T120gross = zeros(NumberWindow,1);
Returnvec3T120gross = zeros(NumberWindow,1);
Returnvec4T120gross = zeros(NumberWindow,1);
Returnvec1T120net = zeros(NumberWindow,1);
Returnvec2T120net = zeros(NumberWindow,1);
Returnvec3T120net = zeros(NumberWindow,1);
Returnvec4T120net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T120(t,:) = onevec./N;
    Returnvec1T120gross(t,1) = wvec1T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T120(t,:)'-wplus,1);
    end 
    Returnvec1T120net(t,1) = (1+Returnvec1T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T120(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T120gross(t,1) = wvec2T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T120(t,:)'-wplus,1);
    end 
    Returnvec2T120net(t,1) = (1+Returnvec2T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T120(t,:) = wghat;
    Returnvec3T120gross(t,1) = wvec3T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T120(t,:)'-wplus,1);
    end 
    Returnvec3T120net(t,1) = (1+Returnvec3T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T120(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T120gross(t,1) = wvec4T120(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T120(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T120(t,:)'-wplus,1);
    end 
    Returnvec4T120net(t,1) = (1+Returnvec4T120gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T120gross = 12*(mean(Returnvec1T120gross)-(gam/2)*var(Returnvec1T120gross));
EU2T120gross = 12*(mean(Returnvec2T120gross)-(gam/2)*var(Returnvec2T120gross));
EU3T120gross = 12*(mean(Returnvec3T120gross)-(gam/2)*var(Returnvec3T120gross));
EU4T120gross = 12*(mean(Returnvec4T120gross)-(gam/2)*var(Returnvec4T120gross));
EU1T120net = 12*(mean(Returnvec1T120net)-(gam/2)*var(Returnvec1T120net));
EU2T120net = 12*(mean(Returnvec2T120net)-(gam/2)*var(Returnvec2T120net));
EU3T120net = 12*(mean(Returnvec3T120net)-(gam/2)*var(Returnvec3T120net));
EU4T120net = 12*(mean(Returnvec4T120net)-(gam/2)*var(Returnvec4T120net));

%%% T = 240
T = 240;
NumberWindow = Ttotal-T;
%Initialize vectors
wvec1T240 = zeros(NumberWindow,N);
wvec2T240 = zeros(NumberWindow,N);
wvec3T240 = zeros(NumberWindow,N);
wvec4T240 = zeros(NumberWindow,N);
Returnvec1T240gross = zeros(NumberWindow,1);
Returnvec2T240gross = zeros(NumberWindow,1);
Returnvec3T240gross = zeros(NumberWindow,1);
Returnvec4T240gross = zeros(NumberWindow,1);
Returnvec1T240net = zeros(NumberWindow,1);
Returnvec2T240net = zeros(NumberWindow,1);
Returnvec3T240net = zeros(NumberWindow,1);
Returnvec4T240net = zeros(NumberWindow,1);
%Rolling window loop
for t = 1:NumberWindow
    %Inputs
    Xis = R(t:t+T-1,:);
    Xoos = R(t+T,:);
    muhat = mean(Xis)';
    sigmahat = cov(Xis,1);
    invsigmahat = inv(sigmahat);
    wghat = invsigmahat*onevec./sum(invsigmahat*onevec);
    %1) EW
    wvec1T240(t,:) = onevec./N;
    Returnvec1T240gross(t,1) = wvec1T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec1T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec1T240(t,:)'-wplus,1);
    end 
    Returnvec1T240net(t,1) = (1+Returnvec1T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %2) EWRF
    wew = onevec./N;
    wvec2T240(t,:) = (wew/gam)*(wew'*muhat)/(wew'*sigmahat*wew);
    Returnvec2T240gross(t,1) = wvec2T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec2T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec2T240(t,:)'-wplus,1);
    end 
    Returnvec2T240net(t,1) = (1+Returnvec2T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %3) SGMV
    wvec3T240(t,:) = wghat;
    Returnvec3T240gross(t,1) = wvec3T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec3T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec3T240(t,:)'-wplus,1);
    end 
    Returnvec3T240net(t,1) = (1+Returnvec3T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
    %4) SMV
    wvec4T240(t,:) = (1/gam)*invsigmahat*muhat;
    Returnvec4T240gross(t,1) = wvec4T240(t,:)*Xoos';
    wplus = zeros(N,1);
    if t == 1
       Turnover = 0;  
    else
        for i = 1:N
            wplus(i,1) = wvec4T240(t-1,i)*(1+R(T+t-1,i));
        end 
        Turnover = norm(wvec4T240(t,:)'-wplus,1);
    end 
    Returnvec4T240net(t,1) = (1+Returnvec4T240gross(t,1))*(1-(cbps/10000)*Turnover)-1;
end
EU1T240gross = 12*(mean(Returnvec1T240gross)-(gam/2)*var(Returnvec1T240gross));
EU2T240gross = 12*(mean(Returnvec2T240gross)-(gam/2)*var(Returnvec2T240gross));
EU3T240gross = 12*(mean(Returnvec3T240gross)-(gam/2)*var(Returnvec3T240gross));
EU4T240gross = 12*(mean(Returnvec4T240gross)-(gam/2)*var(Returnvec4T240gross));
EU1T240net = 12*(mean(Returnvec1T240net)-(gam/2)*var(Returnvec1T240net));
EU2T240net = 12*(mean(Returnvec2T240net)-(gam/2)*var(Returnvec2T240net));
EU3T240net = 12*(mean(Returnvec3T240net)-(gam/2)*var(Returnvec3T240net));
EU4T240net = 12*(mean(Returnvec4T240net)-(gam/2)*var(Returnvec4T240net));

Table4Panel30IND = ...
[EU1T60gross,EU1T60net,EU1T120gross,EU1T120net,EU1T240gross,EU1T240net;...
EU2T60gross,EU2T60net,EU2T120gross,EU2T120net,EU2T240gross,EU2T240net;...
EU3T60gross,EU3T60net,EU3T120gross,EU3T120net,EU3T240gross,EU3T240net;...
EU4T60gross,EU4T60net,EU4T120gross,EU4T120net,EU4T240gross,EU4T240net]

Table4Complete = [Table4Panel10MOM;Table4Panel16ANOM;Table4Panel25SBETA;...
    Table4Panel25SBTM;Table4Panel25OPINV;Table4Panel30IND]